Stochastic extinction of epidemics in large populations and role of vaccinations 
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We investigate stochastic extinction in an epidemic model and the impact of random vaccinations 
in large populations formulated in terms of an optimal escape path. We find that different random 
vaccination strategies can have widely different results in decreasing expected time till extinction, 
for the same total amount of vaccines used. Vaccination strategies are considered in terms of two 
parameters: average frequency of vaccinations, given by 7, and the amplitude of the vaccinations, 
e, where e < 1 refers to the proportion of the population being vaccinated at some particular 
instant. It is found that while the average number of individuals vaccinated per unit time, je, 
is kept constant, the particular values of 7 and e can play a highly significant role in increasing 
the chance of epidemic extinction. The findings suggest that expected time till extinction can be 
significantly shortened if less frequent vaccinations occur in larger groups, corresponding to low 7, 
high e strategy. 
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While deterministic and network models have sug- 
gested certain strategies for epidemic control and pre- 
vention [l], 0], stochastic models are needed to account 
for many of the important features of epidemics, not the 
least of which is stochastic extinction ji ,0, Q • From the 
general theory of finite Markov chains fj| , it was shown 
that in stochastic models the probability of extinction is 
equal to one in the asymptotic time limit. Numerical 
0, El and analytic [3, EH comparisons of stochastic 
and deterministic models have been performed and con- 
firmed that extinction was inevitable in the presence of 
stochastic effects. The numerical results hold for very 
small amplitude noise as well as real finite noise. This 
is in contrast to deterministic SIS or SIR models which 
result in equilibrium endemic presence of infectives [l| 
for an appropriate choice of parameters. It is clear that 
stochastic effects may result in very different dynamics 
from deterministic models, particularly when extinctions 
occur. Since time delivery of vaccines into populations 
is often not a synchronized process, it is important to 
take stochastic effects into account when considering dif- 
ferent vaccination strategies and their impact on extinc- 
tion [3, Q3, E3|- Some numerical comparisons exist on 
the impact of pulsed versus random vaccination strate- 
ies in increasing the probability of stochastic extinction 
However, analytic comparisons of the role of various 
random vaccination strategies in stochastic extinction in 
large populations have not been previously carried out. 

Here we use a deterministic SIS model, with small stan- 
dard deviation Gaussian noise added to model the effect 
of stochastic fluctuations on the dynamics of a large pop- 
ulation. We then consider the effect of different random 
vaccination strategies on the probability of extinction in 
large populations. Our approach differs from the previ- 
ously used Markov models of stochastic epidemics. These 
models have used asymptotic approximations to obtain 
mean extinction times and quasi-stationary distributions 



11, 14, [H, QjiJ . Our method is based on finding the 
optimal escape path that occurs whenever a system ex- 
periences a large and rare stochastic fluctuation from its 



equilibrium state [17j, LL8|, [19j. In this case, the proba- 



bility densities of different trajectories during extinction 
are very different, with highest probability of extinction 
occurring along the optimal escape path. The problem 
of a large fluctuation for stochastic Markov chains was 
treated in [19], where an optimizing action functional is 
used to find the optimal escape path for a particle in a 
well. Here we apply the path integral approach to find the 
optimal escape trajectory that corresponds to stochastic 
extinction of an SIS epidemic when the additive noise is 
sufficiently small. The additive noise, while being less re- 
alistic than multiplicative noise that is sometimes used in 
disease models, allows for a more tractable model, while 
retaining the same qualitative results, with regards to op- 
timal vaccination strategies. In our model, we obtain a 
closed form solution for the stochastic fluctuation needed 
to push the system along the escape trajectory, thereby 
computing the probability of extinction [l3|. The opti- 
mal escape approach simplifies the analysis by allowing 
one to construct a particular deterministic trajectory, or 
optimal escape path, of infectives in the cases that ex- 
tinction occurs. We then use the optimal escape path to 
estimate further the effect that vaccinations have on the 
probability of extinction from a fully endemic state. 

While spatial inhomogeneities may increase stochas- 
tic fluctuations leading to faster extinctions (3|, in gen- 
eral the size of the population is key in determining the 
expected time to extinction For sufficiently large 

populations, normalized fluctuations approach a Gaus- 
sian distribution and scale as one over the square root of 
the population size [20]. We therefore focus on the SIS 
model where the population size, N, is large, so that the 
stochastic fluctuations are relatively small. We model 
these stochastic effects by adding a small standard devi- 
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ation Gaussian noise term, f(t), to the SIS model, 

^=H-f3IS + 8I-iiS + f(t) (1) 

?L=J3IS-8I- ill - f(t) (2) 

where 5 and / are susceptible and infective populations, 
respectively, 8 is a recovery rate, (3 contact rate, and \i is 
the birth/death rate. 

All the variables in Eqns. Jl]) and (J2]) have been 
scaled by the total population size, TV, so that S + I = 
1. For (3 > S + [i, the equilibrium solution, I eq — 
((3 — (8 + /i)) / j3 corresponds to a stable endemic state 
and I ex — is the unstable disease free equilibrium 
(DFE) corresponding to extinction. For smaller values 
of (3, the epidemic dies out even in the absence of any 
stochastic effects, since the DFE is asymptotically stable. 
The term, f(t), denotes random stochastic fluctuations 
between the susceptibles and the infectives, such that 
the total population size is conserved. It is assumed to 
be uncorrelated with zero mean and standard deviation 
D. 

The Langevin formulation, corresponding to Eqns. ([I]) 
and most closely agrees with numerical simulation 
when the disease parameters are near threshold, where 
the equilibrium number of infectives, I eq , is significantly 
smaller than susceptibles. Numerical results, however, 
show good agreement with the Langevin or the Fokker- 
Planck approaches when the equilibrium level of infec- 
tives is as large as one third of the total population (cor- 
responding to A = 1/(1 — I eq ) = 1.5) [5j. The analysis 
in the present paper does not impose any restriction on 
the parameters that determine the equilibrium levels of 
infectives. There is however a constraint on the size of 
stochastic fluctuations, which should be sufficiently small 
for the optimal escape analysis to be applicable. Because 
the size of stochastic fluctuations is determined by the 
size of the total population, the present results should 
also apply for cases that are close to the threshold, pro- 
vided the total population size is sufficiently large. 

Since the stochastic fluctuations are small, the solu- 
tion will spend most of its time close to the equilibrium 
state, I eq - However, given a sufficient amount of time, 
stochastic fluctuations will build up in such a way that 
the number of infectives will go to I ex — 0, leading to 
disease extinction. In large populations, this is a highly 
improbable event, where the length of time to extinction 
has an exponential dependence on the size of the pop- 
ulation [ll|, El. We are interested in how a string of 
vaccinations increases the probability of extinction of an 
epidemic due to stochastic effects. In this case, "vaccina- 
tions" can also stand for various preventative strategies, 
such as wearing face masks, or avoiding all contact with 
the infectives. Suppose at any interval of time, At, there 
is a probability, 7 At, of an occurrence of a vaccination, 



whereby some proportion of the population, e <C 1 is 
vaccinated. The string of vaccinations is then given by 
a sequence of Poisson distributed pulses of amplitude e 
and average number of vaccinations per unit time is "fe. 

Suppose that at some time, t — tj, a number, Nj, of 
individuals are vaccinated, so that Nj/N = e at t = tj, 
and N is the total number of people in the population. 
If each individual acts independently of everybody else, 
than the vaccination process is Poisson distributed with 
e = 1/N. However, it often happens that individuals are 
vaccinated in groups. One example would be when col- 
leges offer vaccinations to all of the students over some 
short period of time, usually involving a few days. In 
this case, different vaccination cites can act as indepen- 
dent agents making a decision on when to offer group 
vaccinations, so that the vaccination process has a Pois- 
son distribution with e = Nj/N, where Nj is the size of 
the group that a particular cite chooses to vaccinate at 
time tj. For simplicity we assume that the size of the 
group being vaccinated is the same at all times (so that 
e is a constant), however the results can be generalized 
to variable group sizes. 

Since vaccinations confer immunity, the number of 
these immune individuals will only change due to death, 
ji. The time-evolution of immune individuals is then 
given by, 

n n 

3=1 3=1 

where ((t) is Poisson distributed, corresponding to a 
string of random vaccinations occurring at various times, 
{ti, t%, ...tn}, with t n < t. Since the total population 
stays constant, I{t) + S(t) + ((t) — 1, we substitute for 
S in Eqn. Q, to get 

^- = al-(31 2 + f(t)-(3iat), (4) 

where a = (3 — 8 — fi, Neglecting the last term, the above 
equation corresponds to a noisy logistic model, used in a 
variety of contexts, such as chemical reactions, transmis- 
sion of rumors, and population growth [l(j]. Equation (f4]) 
also corresponds to the Langevin equation of a particle 
trapped in an over-damped potential well centered at I eq 
with a potential maxima at I ex . 

In the absence of vaccinations, £ = 0, the Gaussian 
noise term, f(t) causes fluctuations about the equilibrium 
state. Given enough time, random noise fluctuations will 
combine in such a way as to drive the particle towards I ex , 
resulting in extinction. For small standard deviation, D, 
of the noise term, /(tLthis extinction will occur along an 
optimal escape path (lR[2lj. For uncorrellated Gaussian 
noise: < f(t)f(tf) >= 5{t — ti), the probability of optimal 
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escape is then given by, 



P[r, 



, SC | (0) = exp[ 



1 

2D 



f opt (t) 2 dt] = exp[-n^/D], 

(5) 

where f op t denotes the stochastic fluctuations that occur 
when the trajectory moves along the optimal escape path 
[l3|, for which H®) is minimized. It can be seen that for 
sufficiently small standard deviation of noise, f opt D, 
the probability of extinction along any other non-optimal 
trajectory becomes negligible. The superscript on the 
P[I esc }^ and ft (0) terms in Eqn. JU indicates that that 
this expression is a solution for the probability of an op- 
timal escape path in the absence of any vaccinations. 

Since any vaccinations should increase the probabil- 
ity of epidemic extinctions, using perturbation theory we 
will obtain an additional correction term, TZ^\ that is 
a direct result of a series of random vaccinations. The 
optimal noise, f op t(t), and the resultant optimal escape 
path is found by minimizing the Gaussian noise over the 
escape trajectory, 



probability of extinction can be obtained as a first order 
correction to TZ (0) in Eqn. © HH, 



f3I esc (t)f opt (t)((t)dt. (9) 



The above equation gives a correction to 1Z that can 
be integrated over the Poisson probability distribution, 
7\ [£(£)], of vaccinations to find the increase in chance 
extinctions. The correction term, in the above equa- 
tion is small compared to TZ^°K However, it can still be 
large compared to the standard deviation of noise, D 
(see Eqn. ([5])), so that vaccinations can significantly in- 
crease chance extinctions. Intuitively, vaccinations in- 
crease the probability of extinction by decreasing the 
amount of stochastic fluctuations needed to push the in- 
fectives along the optimal escape path. 

To find the increase in probability of extinction due to a 
probabilistic sequence of vaccination pulses, we evaluate 
7^ 1 [C] with respect to different realizations of £(i). The 
probability of epidemic extinction along the optimal path 
in the presence of vaccinations is then given by 
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L(i,i)dt. 



(6) 



P[Iesc] =Po I CXpf-K^/D 

J — CO 



7>ac(*)pc(*), (io) 



This is equivalent to minimizing action over a trajectory 
in an Euler-Lagrange system, where the Lagrangian is 

given by L(I,I) = (j - al + j3I 2 ^J . The optimal es- 
cape path then corresponds to the deterministic trajec- 
tory given by the Euler-Lagrange equations, dL/dl = 

jdt. The optimal escape path is a hetero- 



where P = P[Ie SC ] (0 > ■ We are interested in a sequence 
of vaccinations given by a Poisson distribution, with a 
specific realization given by Eqn. |(3]). For a Poisson dis- 
tribution, the probability density, V([C(t)]T>((t), of any 
specific realization of n random vaccinations over the ex- 
tinction interval, 2T, is given by: dt\/2T...dt n /2T, Using 
clinic orbit connecting I eq to I ex , The steady states E q ns ; © and ©- (jjfr, and rescuing time as t -> a 2 t, we 
,I e q and lex, are both saddles in the conservative Euler- 



Lagrangian system given by Eqn. |(6]). The integral 
in Eqns. © and |(6]) is taken for t from — oo to oo, 
which corresponds to the motion along the heteroclinic 
orbit of the system connecting the saddle points. Us- 
ing Euler-Lagrange equations of motion, and conserva- 
tion of energy, we solve for the optimal escape path, 
{I eS c(t) , Iesc(t)} , in the absence of vaccinations. After 
solving for I esc and integrating, we get the trajectory of 
infectives as a function of time for the most probable path 
of extinction, 



get the increase in the chance of extinction, conditional 
on n vaccinations of amplitude e, 



P[Ie 



J f exp I fix J <t>{t)Y^e i 1 



(t-^dt 



dti dt n 
2T "' 2T 

(11) 



where fl = fi/a , T = a T , x is a function of parameters, 



2a 3 ( e@\ _ 2e ( (J3 - 5 - fj,)' 



D 



(12) 



exp (a 2 t) 
1 + exp (a 2 t) 



Using Eqn. © and f op t(t) = I eac (i) - al e 
solve for the optimal noise: 



fo P t{t) 



2a^_ 

T 



exp (a 2 t) 
(1 + exp {a 2 t)f 



(7) 



we 



(8) 



and 4>{t) is the scaled I eS c(t) fopt(i) variable obtained from 
Eqns. and ®, 



exp (2t) 
{l + exp(t)f 



(13) 



Th optimal path is perturbed in the presence of vacci- 
nations, £(£)• F° r small £ (£ <C f op t), however, the per- 
turbation is small and the effect of vaccinations on the 



Since the occurrence of any of the pulses over the interval 
is independent of the other pulses, Eqn. I|lip can be 
rewritten as 

P[Iesc\n] _ 



exp 



)n 
■ 

(14) 
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The above equation gives the probability of escape when 
the number of pulses during the escape interval is n. Us- 
ing a Poisson probability distribution, and summing over 
all possible n, the total escape probability is then given 
by, 



Po 



(15) 



where A n is given by the right hand side of Eqn. (fl4|) . If 
the average number of pulses per unit time is 7, then n = 
27T in the above equation. The sum in Eqn. lfT5|) is an 
expansion of the exponential function, exp [— (1 — A) n]. 
Taking the log on both sides and scaling by a 2 /S, the 



LHS is: 5 



In 



; and the RHS side becomes: 



1 — exp 



jJLX 



<P(t)e 



>dt 



ds (16) 



In the above equation, the limits of integration, ±T have 
been extended to infinity, since the optimal escape tra- 
jectory is along the heteroclinic orbit. Figure Q] plots the 
scaled logarithmic increase in escape probability, S as a 
function of x for p — 0.7, 1, and 2, given by curves (b)- 
(d). As p increases, 2 asymptotes fast to the upper limit 
given by curve (e) in the same figure [12] ■ Letting p — > 00 




Figure 1: Scaled Logarithmic increase in extinction probabil- 
ity, S, as a function of x. (a) Approximation for low x, given 
by Eqn. dHJ) , (b) p = 0.7, (c) p = 1, (d) p = 2, (e) Limit of 
large p, given by Eqn. lfT7|l . (f) Asymptotics valid for large p 
and large x, given by Eqn. (f20|) . 

in Eqn. lfT6|) . the expression corresponding to curve (e) 
in Fig. [l] is, 



(1 — exp [x(/)(s)]) ds 



(17) 



4>(s) is a bounded function with a global maxima at sq = 
In 2, and asymptotically approaching zero as s — * ±00. 
For small x, Eqn. lfl7|) can be approximated as 



<fi(s)ds 



x 

r 







(18) 



This line, 3 = x/2, is plotted in Fig. [H curve (a) and 
provides a good approximation for all values of p, while 



x is sufficiently small. Since x has a linear dependence on 
the vaccination amplitude, e, and E = 7- In I L p " cj J is 
scaled by the average vaccination frequency, 7, we have, 
using Eqns. (O, and flU) , 

fP[Iesc(t)}' 



In 



V Po 



oc 67; 







(19) 



From Eqn. (|19|) . it is clear that increase in extinction 
probability has an exponential dependence on the av- 
erage number of vaccinations per unit time, given by 
67. However, this is only true in the range of smaller 
x, whereby the linear approximation plotted in Fig. [lj 
curve (a) is valid. At higher x, S has a nonlinear depen- 
dence on x (see Fig. [TJ, suggesting that increasing e will 
be more effective in increasing the extinction probability 
than increasing 7. In this range, it is therefore not just 
the average vaccinations per unit time, ej, but the am- 
plitude of these vaccinations, or the degree of fluctuation 
about the average, that is important in decreasing time 
till extinction. 

For large x, we can use Laplace's method [23] to derive 
an asymptotic approximation for Eqn. lfl~7j) . 



2tt 



_ e x<p(to) x 



V -atf"(*o) 
l/(d^M 4 )(i ) 5[0'"(t o )] 2 



x \ 8[<f>"(to)Y 24[0'%)] c 



(20) 



where 4>{to) and its various derivatives are evaluated us- 
ing Eqn. (JT3J) with t = In 2. 

The curve given by Eqn. lf20l) is plotted in Fig. [U 
curve (f). As can be seen from the figure, Eqn. f20|) 
gives an upper limit on the increase in extinction prob- 
ability due to vaccinations. From Eqn. (|20|) . at higher 
x, In (P[I e sc(t)]/Po) has an exponential dependence on 
e and only a linear dependence on 7. It therefore fol- 
lows that in this range of parameters, keeping the average 
number of vaccinations per unit time fixed, the strategy 
of delivering high amplitude lower frequency vaccinations 
is much more effective in increasing the probability of 
stochastic extinctions in a fully blown epidemic. This 
perhaps makes sense in the context of stochastic extinc- 
tions, when one considers that in the absence of any vac- 
cinations, the time till extinction depends on the size of 
stochastic fluctuations relative to the size of the popula- 
tion. Random vaccinations can themselves be considered 
as a source of positive stochastic fluctuations, adding to 
the stochastic fluctuations in disease transmission. Since 
less frequent, higher amplitude vaccinations correspond 
to a greater standard deviation of vaccinations, they add 
more to the stochastic fluctuations of the whole popula- 
tion. This results in lesser expected time till extinction, 
compared to other random vaccination strategies that use 
the same number of vaccines. 
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The analysis in this paper suggests that depending on 
parameters (given by birth, recovery and contact rates), 
it can be far more effective to vaccinate individuals in 
groups, rather than allowing each individual to make an 
isolated decision. This can be understood as follows: if 
each person vaccinates independently of everybody else 
than e takes the smallest possible value of 1/N, corre- 
sponding to the amplitude of an individual vaccination. 
However, if a group of size Nj vaccinates at approxi- 
mately the same time, then the amplitude, e is given by 
Nj/N. It follows that even if or the average num- 
ber of individuals vaccinating per unit time is the same 
in both cases, the exptected time till extinction may be 
significantly shortened. 
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